Summary

Adds conservative regridding for grids that aren't 1D-separable — a capability the existing .conservative method cannot express:

• Curvilinear grids — accepts 2D lat[i, j] / lon[i, j] coordinate variables (ORCA, tripolar, rotated regional).
• Unstructured meshes — ConservativeRegridder.from_polygons(...) accepts arbitrary shapely polygon arrays (MPAS, ICON, finite-element).
• Grid → polygon aggregation — regrid structured data onto country shapes or watershed boundaries, same from_polygons path. Auxiliary target coordinates (region_id, labels, scalar metadata) ride along with the output.
• Correct spherical areas on lat/lon — opt-in spherical=True kwarg applies an analytic Lambert cylindrical equal-area projection of cell edges before intersecting. Matches the factored path's sin-weighted accuracy at the same fast-path cost.
• Antimeridian support — structured lat/lon grids whose longitude spans the dateline are unwrapped automatically. For polygon input, from_polygons(..., periodic=True) unwraps polygon rings that cross ±180°.

Internals:

• ConservativeRegridder class — reusable regridder; builds the sparse area-intersection matrix once, applies to many fields. .T / .transpose() returns the backward regridder sharing cached state.
• .to_netcdf / .from_netcdf — persist the weight matrix plus reproducibility metadata (dim names, shapes, grid ranges, spherical flag, xarray-regrid version, timestamp, schema version) in root attributes, with source / target coord-only groups alongside. Schema-versioned for forward compat.
• Accessor — exposed via .regrid.conservative_2d(...) on xarray DataArrays / Datasets.
• Optional extra — pip install xarray-regrid[conservative-2d] pulls shapely>=2.0 + h5netcdf. Existing users unaffected.

Motivation

.conservative uses Stephan Hoyer's axis-factored 1D overlap — fast and elegant but strictly rectilinear (1D-separable). Users with curvilinear ocean models, unstructured climate meshes, or grid-to-region aggregation currently reach for xESMF, which:

• Requires ESMPy (C/Fortran via conda-forge only; no Windows support)
• Doesn't support dask.distributed
• Is a heavy dep cascade

.conservative_2d fills the gap with pure shapely + sparse. Works on Windows, composes with distributed dask, installs with a single pip install shapely.

Naming follows xESMF's short-method-name convention. "_2d" distinguishes full 2D cell-polygon intersection from the 1D axis-factored approach and doesn't prejudge whether the grid is curvilinear or unstructured.

Design notes

• Rectilinear fast path: when both grids are axis-aligned rectangles, skip GEOS clipping entirely and compute intersection areas analytically from the bounds. ~30× faster than the GEOS path at 1M cells.
• STRtree candidate filter uses bounding-box overlap, not the GEOS intersects predicate, by default — much cheaper and the subsequent area > 0 filter drops any false positives. from_polygons(..., predicate_filter=True) restores the predicate filter for pathological inputs (thin diagonals with loose bboxes).
• Threaded GEOS: when curvilinear, shapely 2's GIL-free shapely.intersection over candidate pairs is parallelized across a ThreadPoolExecutor. ~3-4× on 8 cores.
• Sparse weight matrix: stored as sparse.COO; apply goes through sparse.matmul inside xr.apply_ufunc(dask=\"parallelized\"). Cached pre-sorted transpose avoids a per-call re-sort.
• NaN semantics match existing .conservative: two-pass value + mask matmul, nan_threshold with the same interpretation.
• Coverage mask applied regardless of skipna so uncovered target cells (domain boundaries, polygon holes) always produce NaN.
• Dtype preserved: float32 in → float32 out. sparse matmul promotes to float64 internally; cast back in _apply_core.

Test plan

• Structured rectilinear: conservative_2d exactly reproduces .conservative (planar) to 1e-12.
• Dask + chunked input: identical results to non-dask; spatial chunks rechunked automatically.
• Curvilinear target (rotated 2D): mass conservation at 5.24e-16 rel.
• Unstructured mesh: 50k-cell Voronoi round-trip at 1e-15 rel.
• NaN handling: strict vs lenient thresholds produce the expected ordering.
• .T transpose on rectilinear and curvilinear; constant field round-trip exact on aligned grids.
• Spherical mode: matches factored path on lat/lon grids to within grid quadrature.
• Antimeridian: constant field through a dateline-crossing structured grid is preserved exactly; from_polygons(periodic=True) correctly routes samples across the seam.
• Auxiliary target coords (non-dim + scalar) attached on from_polygons output.
• netCDF save/load: bit-identical output on reload; schema mismatch rejected.
• Dtype preservation: float32/64 roundtrip as expected.
• Empty target / polygon holes → NaN (regression-guarded).
• 33-test test_conservative_2d.py module passing on this branch.

Docs

Three demo notebooks under docs/notebooks/demos/:

• demo_conservative_2d_curvilinear.ipynb — rotated 2D target.
• demo_conservative_2d_regions.ipynb — structured-grid → region-polygon aggregation.
• demo_conservative_2d_unstructured.ipynb — Voronoi mesh via from_polygons + .to_netcdf / .from_netcdf.

Notebooks are output-stripped — reviewers / docs build execute them fresh.

Follow-ups (out of scope for this PR)

• Great-circle (s2) geometry: WIP on spherely-integration branch (PR WIP: Add s2 manifold to conservative_2d (stacks on #70) #71) — drops in manifold=\"s2\" backed by the optional spherely dependency.
• Grid-to-region aggregation (like xagg): natural extension of from_polygons but wants a dedicated API surface for loading region polygons.

🤖 Generated with Claude Code